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INVESTIGATION  OF  MOLECULE-SURFACE  INTERACTIONS 
WITH  OVERTONE  ABSORPTION  SPECTROSCOPY  AND 
COMPUTATIONAL  METHODS 


1 .  INTRODUCTION 

1 . 1  Objective 

The  energetics  and  dynamics  of  molecule-surface  interactions  play  an  integral 
part  to  chemical  processes  and  chemical  adsorption.  Chemical  reactions  on  the  surface  of  a 
catalyst  arc  a  good  example  of  these  interactions.1  A  specific  set  of  applications  of  molecule- 
surface  interactions  that  is  of  interest  to  the  warfighter  is  in  chemical-biological  defense.  For 
example,  chemical  personal  protective  equipment  relics  heavily  on  a  material’s  ability  to  surface 
bind  harmful  compounds.  Another  type  of  study  of  interest  to  the  Army’s  mission  is  the 
interaction  of  energetic  material  residues  with  various  environmental  surfaces.  These 
interactions  can  have  an  impact  on  detectability  of  the  energetic  materials.  A  number  of  surface 
techniques  have  been  applied  to  this  problem,  such  as  diffuse  reflectance  IR  spectroscopy'  and 
step  scan  FTIR  spectroscopy.1,3  These  techniques  have  been  effective  at  identifying 
decomposition  intermediates  and  products  even  on  the  timescale  of  hundreds  of  nanoseconds. 
Although  these  techniques  are  proven  tools  for  the  study  of  the  chemical  dynamics  of 
decomposition,  details  of  the  molecule/surface  interaction  arc  still  unknown.  Vibrational 
overtones  arc  uniquely  sensitive  to  a  molecule's  local  environment,  and  are  potentially  an 
effective  tool  in  examining  molecule/surface  interactions.  It  is  the  objective  of  this  work  to 
determine  the  most  effective  and  general  theoretical  model  of  molecular-surface  interactions. 
This  will  be  done  using  vibrational  overtone  spectroscopy  in  tandem  with  molecular  mechanics 
and  ab  initio  calculations  to  investigate  molecule-surface  binding.  Due  to  the  importance  of 
1,3,5-trinitrotoluene  (TNT)  as  a  compound  and  the  catalytic  and  optical  properties  of  alumina, 
the  systems  of  neat  TNT  and  TNT  on  an  alumina  surface  serves  as  the  test  cases  for  this  study. 

1 .2  Background 

The  interaction  of  a  molecule  with  its  local  environment  is  reflected  in  its 
vibrational  spectrum,  so  that  the  vibrational  spectrum  can  be  used  as  an  indicator  of  the  local 
environment.  Vibrational  overtone  spectroscopy  is  ideal  for  studying  molecule-surface 
interactions  because  it  involves  higher  vibrational  energy  levels.  When  the  vibrational  potential 
energy  surface  of  a  molecule  is  perturbed  by  the  presence  of  a  surface  or  solvent,  the  higher 
energy  levels  experience  a  more  drastic  change  than  the  lower  levels.  This  makes  sense  in  view 
of  the  fact  that  at  higher  vibrational  energy  levels,  the  atomic  displacement  along  the  vibrational 
coordinate  is  greater,  so  that  the  molecule-environment  interactions  will  be  greater  than  at  lower 
vibrational  energies.  These  levels  are  not  typically  accessed  with  fundamental  absorptions 
because  almost  all  molecules  are  in  the  ground  vibrational  state  at  room  temperature,  so  that 
optical  absorption  only  promotes  molecules  to  the  first  vibrationally  excited  state.  This  energy  is 
quickly  lost  to  equilibrium  with  the  environment  so  that  multiphoton  processes  are  unlikely.  As 
a  result,  the  higher  vibrational  states  are  only  accessible  by  overtone  transitions,  that  is, 
transitions  involving  more  than  one  quantum  of  vibrational  energy.  Although  molecules  are 


7 


typically  modeled  as  harmonic  oscillators,  real  molecules  arc  anharmonic  oscillators  where 
transitions  of  more  than  one  vibrational  quantum  are  possible/  As  a  result,  the  anharmonic 
description  of  molecular  vibrations  is  necessary  for  this  study. 

Most  theoretical  models  use  the  harmonic  description  of  molecular  vibrational 
because  it  is  computationally  easier,  and  because  for  studies  of  fundamental  transitions  it  is 
normally  adequate.  However,  there  are  additional  differences  between  the  anharmonic  and 
harmonic  description  of  molecular  vibrations  that  show  that  the  harmonic  description  is  not 
sufficient  for  this  study.  In  the  harmonic  description,  the  potential  energy  surface  goes  to  infinity 
at  large  atomic  separation  and  does  not  account  for  the  breakage  of  bonds.  The  result  is  that  the 
calculated  vibrational  frequencies  for  a  harmonic  oscillator  arc  integer  multiples  of  the 
fundamental  frequency  that  go  up  infinitely.  In  real  anharmonic  systems,  the  separation  between 
energy  levels  tends  to  zero  as  the  vibrational  energy  approaches  the  dissociation  limit.  This 
means  that  the  overtone  transition  frequency  will  not  be  an  integer  multiple  of  fundamental 
transition  frequency,  and  that  a  perturbation  is  likely  to  increase  the  deviation  from  the  frequency 
predicted  by  the  harmonic  treatment.  In  other  words,  the  harmonic  oscillator  description  will 
consistently  overestimate  molecular  fundamental  and  overtone  vibrational  frequencies.  The 
anharmonic  description  is  essential  for  adequate  handling  overtone  vibrations. 

There  are  several  available  approaches  for  modeling  the  interaction  of  a  molecule 
with  a  surface;  however,  only  some  of  these  methods  permit  prediction  of  overtone  frequencies. 
The  simplest  is  molecular  mechanics,  which  operates  on  a  quasi-classical  level.6'9  The  details  of 
the  electronic  states  are  neglected  and  the  energy  of  the  system  is  treated  as  a  sum  of  interactions 
through  chemical  bonds,  such  as  bond  stretches  and  bends,  and  non-bonding  interactions  such  as 
electrostatic  interactions  and  van  der  Waals  energy.  The  overall  description  of  the  molecule  as  a 
sum  of  classical  interactions  is  also  known  as  a  “force  field”.  The  strength  of  this  approach  is 
that  it  is  not  computationally  expensive  to  perform,  so  that  large  systems  can  be  modeled.  A 
major  assumption  of  this  approach  is  that  the  interaction  between  two  atoms  in  one  molecule  will 
be  like  that  in  another  molecule.  This  assumption  leads  to  two  weaknesses  of  this  approach;  (1) 
empirical  data  is  often  necessary  to  build  a  complete  force  field,  and  (2)  accuracy  of  the  results 
depends  on  the  similarity  of  the  system  under  study  to  the  systems  that  provided  the  empirical 
data.  Predicted  properties,  such  as  vibrational  spectra,  are  only  accurate  if  the  source  of  the 
empirical  data  is  very  much  alike  to  the  model  system.  Although  typical  force  fields  permit  use 
of  the  Morse  Potential  for  describing  bond  stretches  and  anharmonic  vibrations,  the  accuracy  is 
most  likely  insufficient  to  make  reliable  predictions  of  overtone  frequency  shifts  the  can  be  used 
to  infer  the  molecule/environment  interaction. 

The  system  can  also  be  handled  with  ab  initio  quantum  mechanical  or  density 
functional  theory  (DFT)  methods  that  do  not  require  any  empirical  data.  For  ab  initio  methods, 
the  electronic  wavefunction  of  the  system  is  calculated,  and  the  vibrational  potential  energy 
surface  is  directly  calculated  from  the  electronic  wavefunction.  10  For  density  functional 
methods,  the  square  of  the  wavefunction,  that  is,  the  probability  density,  is  calculated  and  used  to 
get  the  vibrational  potential  energy  surface.  For  some  quantum  mechanical  methods  such  as 
Hartree  Fock  (HF),  DFT,  second  and  third  order  Moeller/Plesset  (MP2/MP3),  it  is  possible  to 
obtain  analytic  second  derivatives  and  numerical  third  derivatives  of  the  vibrational  potential 
energy  surface.  This  allows  more  accurate  anharmonic  treatment  of  the  potential  energy  surface. 
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and  calculation  of  the  overtone  frequencies.  However,  these  approaches  are  computationally 
expensive,  so  that  there  is  a  strong  limitation  on  the  size  of  the  theoretical  model  that  can  be 
handled.  Fortunately,  a  number  of  theoretical  modeling  packages  include  the  capability  to 
perform  hybrid  calculations  where  part  of  a  system  is  treated  with  quantum  mechanics  and 
another  part  is  treated  with  molecular  mechanics.  Proper  application  of  the  hybrid  treatment 
permits  accuracy  approaching  that  of  a  purely  quantum  mechanical  method  at  a  fraction  of  the 
computational  cost. 

The  overtone  frequency  serves  as  the  bridge  between  the  molecule-surface 
interaction  model  and  experiment.  A  modeling  approach  will  be  considered  valid  if  the 
predicted  overtone  frequencies  accurately  reproduce  those  measured  from  experiment.  A 
number  of  important  details  in  the  surface/moleeule  interaction  play  a  part  in  the  simulation  that 
w  ill  impact  predicted  overtone  frequencies.  For  example,  the  orientation  of  the  molecule  relative 
to  the  surfaee  will  affeet  how  its  force  Field  is  perturbed,  so  that  each  possible  orientation  will 
have  its  own  spcetral  shift.  If  no  particular  orientation  is  favored,  then  the  observed  shift  in  the 
overtone  spectrum  will  be  a  convolution  of  transition  energies  of  all  possible  orientations.  If  a 
particular  orientation  is  favored,  as  in  a  catalytic  surfaee,  then  the  shift  in  the  overtone  spectrum 
will  reflect  this.  Another  important  detail  is  the  magnitude  of  the  moleeule-surface  binding 
energy.  The  binding  energy  will  govern  transport  dynamics.  The  binding  energy  is  easily 
determined  in  the  model  by  the  difference  between  the  energy  of  the  molecule  and  the  surface 
separated  by  infinite  distance  and  the  energy  of  the  moleeule-surfaee  system  at  the  optimized 
geometry.  It  is  expeeted  that  a  stronger  binding  energy  will  correspondingly  perturb  the 
molecular  vibrations  and  result  in  more  drastie  shifts  in  the  overtone  absorption  frequency 
relative  to  the  isolated  molecule. 


2.  EXPERIMENT  AND  SIMULATION 

2.1  Theoretical  Considerations 

2.1.1  Computational  Approach 

Neat  TNT  and  TNT  molecules  on  an  alumina  surface  are  modeled  utilizing 
Gaussian  2003  version  B.04  on  a  Sunstation  V,  as  well  as  a  LAN  Parallel  cluster  utilizing 
48  Intel  processors,  and  Gaussian  2009,  version  A. 02.  The  Sunstation  system  has  32  Gbytes  of 
RAM  memory  and  eight  processors;  however,  this  implementation  of  Gaussian  2003  does  not 
include  the  LAN  parallel  processing  eode.  Each  node  in  the  Gaussian  2009  cluster  has  eight 
processors  and  eight  Gbytes  of  RAM.  To  calculate  the  overtone  absorption  frequencies  and 
intensities  of  TNT  interacting  with  a  surfaee,  anhannonie  frequency  analysis  and  the 
Harmonically  Coupled  Anhannonie  Oscillator  (HCAO)  model  are  applied.  In  the  HCAO  model, 
eaeh  C-H  unit  within  the  moleeule  is  treated  as  a  diatomic  Morse  oscillator  that  is  harmonically 
coupled  to  the  other  C-H  units. 
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The  first  task  is  to  determine  an  appropriate  level  of  theory  for  performing 
anharmonic  frequency  and  HCAO  analysis  on  the  TNT  molecule.  The  criteria  for  selection  are 
minimal  computational  expense  and  accuracy.  Accuracy  is  determined  by  comparison  of  the 
calculated  molecular  dimensions  and  fundamental  vibrational  frequencies  to  published  X-ray 
data  as  well  as  spectroscopic  data.  Anharmonic  frequency  calculations  are  run  using 
Hartree/Fock  (HF),  second  order  Moeller/Plesset  (MP2),  and  density  functional  theory  (DFT), 
methods  that  have  available  analytic  second  derivatives  of  the  potential  energy  surface.  Basis 
sets  up  to  6-31 1+G(d,p)  arc  evaluated.  More  complex  basis  sets  that  have  more  diffuse  functions 
are  not  attempted  for  anharmonic  frequency  calculations  due  to  computational  limitations. 
Attention  is  limited  to  the  TNT  C-H  stretches  because  they  occur  at  higher  frequencies  and  the 
resulting  overtone  absorptions  are  accessible  from  wavelengths  produced  by  the  optical 
parametric  oscillator  (OPO)  laser  system.  For  HCAO  analysis,  DFT  with  the  B3LYP  functional 
is  used  as  a  starting  point  based  on  the  inaccuracy  of  frequency  calculations  with  HF  and  the 
computational  cost  of  utilizing  higher  levels  of  theory  such  as  MP2. 

The  second  task  is  the  calculation  of  absorption  frequencies  as  a  function  of  local 
environment.  Several  systems  are  investigated  with  the  standard  frequency  analysis  within  the 
Gaussian  code  to  calculate  the  optimal  geometry  followed  by  the  normal  frequencies.  Normal 
mode  analysis  is  then  followed  by  determination  of  anharmonic  corrections.  For  Gaussian,  the 
second  order  force  constants  arc  determined  analytically  to  give  the  cubic  force  constants  that  are 
used  to  determine  the  anharmonic  corrections.  Within  Gaussian  2003,  specifying  anharmonic 
analysis  limits  the  size  of  the  system  that  can  be  studied  because  it  attempts  to  perform  the 
analysis  for  every  mode.  Gaussian  2009  permits  the  user  to  limit  the  scope  of  the  anharmonic 
analysis.  Thus,  results  are  obtained  from  some  hybrid  quantum  mechanical/molccular 
mechanical  models.  Several  systems  arc  studied  with  the  method,  including  single  molecules, 
and  clusters,  where  a  single  TNT  molecule  is  treated  with  quantum  mechanics,  and  the  rest  of  the 
cluster  is  treated  with  the  Universal  Force  Field  (UFF),  allowing  van  der  Waals  and  electrostatic 
interactions  to  perturb  the  quantum  mechanical  portion  of  the  model.  The  electrostatic 
interactions  are  treated  with  (1/R)  dependence,  and  the  van  der  Waals  interactions  are  treated 
with  a  “soft”  cutoff. 

2.1.2  Molecule/Surface  Modeling 

Several  different  approaches  are  evaluated  for  modeling  a  surface,  such  as 
modeling  with  periodic  boundary  conditions  (PBC),  and  hybrid  molecular  mcchanics/quantum 
mechanics  methods.  A  calculation  using  the  PBC  utility  is  attempted  to  determine  the  overtone 
absorptions  in  solid  neat  TNT.  X-ray  crystallographic  measurements  of  the  neat  TNT  crystal 
show  an  eight  member  unit  cell  and  the  spacing  between  the  individual  molecules.12  The 
published  geometry  is  used  as  a  starting  point  for  the  calculation  w'ith  HF  and  the  STO-3G  basis 
set,  which  is  the  lowest  level  method  that  can  be  invoked  by  the  PBC  convention.  Although  half 
of  the  Sunstation’s  resources  are  used  for  this  calculation  (four  processors  and  16  Gbytes  of 
RAM  memory),  the  calculation  is  terminated  after  20  days  without  completing  the  geometry 
optimization.  This  approach  is  not  feasible  with  the  computational  resources  at  hand. 
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The  next  attempt  at  modeling  a  cluster  of  TNT  molecules  uses  the  AM  1 
Hamiltonian,  which  is  a  semi-empirical  method.  The  semi-empirical  methods  arc  much  less 
computationally  expensive  than  ab  initio  methods,  and  complete  on  a  Windows  desktop  machine 
for  clusters  with  three  and  eight  members.  The  results  of  the  geometric  optimization  did  not 
compare  well  with  published  X-Ray  data,  and  the  difference  between  the  resulting  structure  and 
the  published  TNT  crystal  structural  data  suggests  that  a  very  large  cluster  is  necessary  to 
replicate  the  published  experimental  data.  The  next  approach  is  to  speed  up  the  PBC  calculation 
by  providing  a  unit  cell  optimized  with  a  semi-empirical  method.  This  approach  is  not  workable 
because  the  results  show  significant  deviation  from  the  X-ray  data. 

Calculation  using  the  hybrid  quantum  mcchanics/molecular  mechanics  (QM:MM) 
utility  of  Gaussian  is  also  attempted  using  combinations  of  HF.TJFF  and  B3LYP:UFF.  To  avoid 
application  of  the  anharmonic  analysis  to  atoms  within  the  model  of  the  surface,  Gaussian  2009 
is  used.  Within  the  newer  version  of  the  code,  the  whole  system  geometry'  is  optimized  with 
UFF.  and  then  all  atoms  except  the  TNT  molecule  are  excluded  from  further,  higher  level 
geometry  optimization  using  HF  or  B3LYP.  The  TNT  molecule  alone  undergoes  frequency 
analysis;  whereas,  the  adjacent  surface  model  serves  as  a  perturbation. 

2.1.3  Theory  of  the  HCAO  Model  and  Calculation  of  C-H  Overtone  Frequencies 

and  Intensities 

We  also  use  the  Harmonically  Coupled  Anhannonie  Oscillator  (HCAO)  model  to 
simulate  overtone  transitions  from  the  target  molecule,  1 ,3,5-trinitro-toluene  (TNT)  interacting 
w  ith  a  surface.  The  HCAO  model  allows  calculation  of  absorption  frequencies  and  intensities  of 
the  overtone  transitions,  and  has  been  successful  in  the  literature  at  describing  the  absorptions  of 
several  molecules,  including  chemical  agents.13"22 

A  key  assumption  of  the  HCAO  model  is  that  for  overtones,  especially  the  Au>2, 
the  X-H  bond  stretching  coordinates  become  almost  completely  independent  of  the  other  internal 
coordinates  of  the  molecule,  c.g.,  bends,  torsions,  etc.  For  the  purely  local  mode  model, 
coupling  between  other  internal  motion  coordinates  is  completely  neglected;  whereas.  HCAO 
allows  for  coupling  between  X-H  oscillators  within  the  molecule,  given  that  they  have  similar 
fundamental  frequencies.  For  TNT,  there  arc  five  C-H  oscillators:  three  oscillators  on  the 
methyl  group,  and  two  oscillators  that  are  part  of  the  aromatic  ring.  The  three  methyl  oscillators 
are  treated  as  kineticallv  and  potentially  coupled,  and  the  two  aryl  oscillators  arc  treated  as 
potentially  coupled  with  the  kinetic  energy  coupling  neglected  due  to  the  physical  separation  of 
the  oscillators. 


Within  the  local  mode  model,  the  C-H  oscillators  within  TNT  are  treated  ns 
diatomic  anhannonie  oscillators.  The  Morse  Potential  approximates  the  potential  energy  surface 
of  a  diatomic  molecule  quite  well,  even  taking  into  account  the  dissociation  limit,  where  the 
Harmonic  Oscillator  model  docs  not.  The  potential  energy  V  as  a  function  of  the  displacement 
coordinate  q  around  the  equilibrium  bond  length,  is 


V{q)  =  Dt,(l  —  exp(-/ft/))2 


(1) 


where  q=r-ro,  the  bond  length  r  minus  the  equilibrium  bond  length  r<>,  De  is  the  dissociation  limit 


measured  from  the  bottom  of  the  potential  energy  surface,  and  p  = 


where  c  is  the 


speed  of  light,  p  is  the  reduced  mass,  cocxe  is  the  anharmonicity  constant,  and  h  is  Planck’s 
constant,  6.626  x  10'34  J  s.  Conveniently,  there  exist  analytic  solutions  to  the  nuclear 
Sehroedinger  for  the  Morse  oscillator.  The  Hamiltonian  of  the  resulting  system  for  the  pure  local 
mode  model  can  be  expressed  as23, 


he 


=  U6)e  |t>)  -  6)eXe(l>2  +  L>)|t>) 


(2) 


where  |u>  denote  the  Morse  Oscillator  wave  functions.  Insertion  of  the  appropriate  vibrational 
quantum  numbers  and  Morse  Oscillator  parameters  give  the  transition  frequencies  in  cm  .  For 
the  coupled  oscillators  in  the  TNT  molecule,  the  Hamiltonian  is  written  with  additional  coupling 
terms,  and  the  nuclear  wavefunction  for  the  coupled  oscillators  are  written  as  products  of  the 
individual  wave  functions  along  each  stretching  coordinate.  The  coupling  correction  is 
considered  as  pairwise  interactions  between  the  absorbing  C-H  oscillator  and  the  rest  of  the 
oscillators  in  the  methyl  group.  The  result  is  Equation  3  for  the  coupled  methyl  C-H  stretches, 


he 


=  «*»,uMI0),|0)*  -«VWy2  +  wMI°)jl0)*  -••• 


•••-Z  <p'Ma)  +«xRl°>/l°)1  -• 

i*j 

•••-Zp'<*(« +«x]k),l0),l0)* 


(3) 


where  u  is  the  vibrational  quantum  number,  and  (0c(i)  denotes  the  Morse  oscillator  harmonic 
frequency  of  the  i‘‘  oscillator,  the  optical  radiation  absorbing  oscillator,  <}>„’  is  the  oscillator 
coupling  parameter  between  the  ith  and  jth  coupled  oscillators,  and  a'  and  a  arc  the  harmonic 
oscillator  raising  and  lower  operators,  where  it  is  assumed  they  act  on  the  Morse  oscillator  wave 
functions  in  the  same  way  as  the  Harmonic  wave  functions.  The  combination  of  raising  and 
lower  operators  extracts  the  energy  differences  between  coupled  oscillators.  Mathematically,  the 

operators  behave  as,  at  |u).  =  VZ>|u-l)  and  a*\v)  =  i>|l>  +  1)  .  For  the  aryl  C-H  stretches,  the 

Hamiltonian  is  the  same,  except  that  the  nuclear  wave  functions  are  products  of  only  two  internal 
coordinate  wave  functions,  and  there  is  only  one  coupling  correction  term.  The  coupling 
correction  term  <t>jj’  can  be  expressed  as, 

%'=  (r# (4) 

where  y,,  is  the  kinetic  energy  coupling,  (Jig  is  the  potential  energy  coupling,  and  coe(i)  and  coc( j >  are 
the  Morse  frequencies  of  the  i,h  and  jlh  oscillator.  The  coupling  constants  arc  related  to  the 
Wilson  G  and  F  matrix  elements  by,19’20,23,24 
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(5) 


y„  = 


i  g, 


(b) 


where  Gjj  and  Fjj  are  the  off-diagonal  matrix  elements.  For  the  coupled  methyl  C-H  oscillators 
located  on  the  same  methyl  group,  yi,  can  he  written  as,  from  the  tabulated  G  elements,25'2 


7;,=- -cos(^) 
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(7) 


where  0  is  the  angle  between  the  methyl  C-H  bonds,  and  me  and  dim  are  the  carbon  and  hydrogen 
atomic  masses,  respectively.  The  calculation  of  the  fundamental  and  overtone  vibrational 
frequencies  is  straightforward.  To  simplify  the  determination  of  force  constants  from  the 
computational  output,  the  potential  energy  surface  is  treated  as  a  Taylor  series  expansion  up  to 
fourth  order: 1 


V(q)  = 


i  c-v 


cq- 
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i  d'y 
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where  C  is  a  constant,  and  the  derivatives  are  the  square,  cubic,  and  quartic  force  constants, 
evaluated  at  the  equilibrium  position,  where  the  linear  derivative  is  0.  The  force  constants  are 
related  to  the  Morse  parameters  coc  and  cocxe  by. 
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where  G,,  is  one  of  the  diagonal  Wilson  G  matrix  elements.  In  the  case  of  a  bond  stretching 
internal  coordinate,  the  G  element  is  0  =  v/m  \  +  1  where  mx  is  the  atomic  mass  of  X 

and  mu  is  the  atomic  mass  of  Hydrogen  in  the  X-H  oscillator.  The  F  and  G  matrices  relate  the 
internal  coordinates  of  the  molecule  to  the  normal  modes. 


To  extract  force  constants  from  the  computational  output,  the  geometry  of  the 
target  molecule  and  its  environment  are  energy  optimized  with  Gaussian.  Then,  the  potential 
energy  is  calculated  for  a  number  of  grid  points  around  the  equilibrium  position.  In  this  study, 
the  bond  lengths  are  varied  +0.30  A  in  steps  of  0.05  A.  for  a  square  grid  pattern.  The  potential 
energy  scans  are  performed  along  two  stretch  coordinates  simultaneously  to  be  able  to  capture 
any  interaction  between  two  coordinates.  If  two  local  modes  are  completely  uncoupled,  then 
motion  along  one  coordinate  has  no  effect  on  the  other.  For  the  methyl  C-H  oscillators,  each 
oscillator  is  scanned  pairwise,  so  that  interaction  between  all  three  C-H  oscillators  is  captured. 
The  aryl  C-H  coordinates  and  an  aryl  and  methyl  C-H  coordinates  are  scanned  together  to  check 
for  coupling.  For  determination  of  the  cubic  and  square  force  constants,  the  potential  energy  as  a 
function  of  motion  along  one  coordinate  is  fitted  to  equation  (8).  and  all  other  coordinates  are 
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held  at  the  equilibrium  position.  The  coefficients  of  the  fit  generated  by  Matlab  R2006  are 
substituted  into  equations  (9)  and  (10)  to  yield  the  Morse  oscillator  parameters.  The  F  matrix 
elements  are  determined  in  a  similar  way,  except  that  the  potential  energy  surface  is  fitted  along 
two  stretch  coordinates  simultaneously.  Approximating  the  coupling  between  oscillators  as 
harmonic,  only  linear  and  square  terms  are  used  for  the  fit,  yielding  the  F  matrix  elements  that 
are  used  for  calculating  the  coupling  corrections  in  equations  (6)  and  (4).  Sufficient  information 
is  then  on  hand  for  determining  the  full  HCAO  predicted  oscillator  frequency. 


To  calculate  predicted  intensities  using  the  HCAO  model,  the  matrix  elements  of 
the  electric  dipole  moment  operator  |<v|/i|p|v|/o>|  need  to  be  calculated  for  each  oscillator  . 
Because  the  dipole  moment  function  p  is  not  known  analytically,  it  is  approximated  with  a 
Taylor  expansion  up  to  the  fifth  order  as  follows  for  the  three  coupled  methyl  C-H  oscillators 
(11).  In  the  case  of  the  aryl  oscillators,  their  dipole  moment  function  is  treated  as  an  expansion 
along  two  stretching  coordinates  with  terms  taken  up  to  the  fourth  order  (12): 
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It  must  be  pointed  out  that  the  higher  order  terms  contribute  significantly  to  the 
intensity  calculation.  To  obtain  the  dipole  moment  derivatives,  the  dipole  moment  magnitude  is 
determined  on  the  same  grid  points  as  the  molecular  potential  energy  in  the  DFT  or  ab  initio 
calculations.  The  derivatives  are  found  in  the  same  way  as  the  force  constants,  except  the  dipole 
moment  points  arc  fitted  to  polynomials  as  denoted  by  equations  (11)  and  (12).  The  resulting 
dipole  moment  function  is  used  for  numerical  integration  and  determination  of  the  transition 
electric  dipole  moment.  For  the  calculation  of  the  transition  dipole  moment,  the  Morse 
Oscillator  wave  functions  arc  used  as  a  basis  set.  For  the  Morse  wave  functions,2*  the  phase 
factor  is  neglected,  and  the  wave  functions  are  written  as  (13), 


V'(<l)u  =  | L*) =  \ ^  expf _  4 exp(- pq )  jA  exp(-  fiq )  2  cxp(-/3q)-(k  - 2v- 1)) 


(13) 


where  k  is  the  ratio  of  coe  to  coexc,  Lv  is  the  Lagucrrc  polynomial  of  the  terms  within  the 
parentheses.  The  harmonic  basis  wave  functions  along  individual  coordinates  arc  written  as. 
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where  ti  is  h/2jt.5  The  wavefunetion  then  used  for  the  overtone  transitions  is  fomied  from  a 
linear  combination  of  the  harmonic  basis  functions.  For  cases  where  the  k  factor  becomes  large 
(k>  175),  that  is,  the  ratio  of  the  harmonic  frequency  to  the  anhamionicity  constant,  the 
exponential  factor  in  the  Morse  wavefunetion  expression  becomes  too  small  to  represent  with 
double  precision  floating  point  numbers.  In  this  case,  the  linear  combinations  of  the  harmonic 
oscillator  wave  functions  are  used  in  the  calculation  of  intensities. 

To  determine  the  relative  transition  intensities,  the  transition  dipole  moment  |pi0|2 
is  calculated  from  the  wavefunetion  expressions  for  the  ground  and  v=2  excited  states  and  the 
dipole  moment  function.  For  the  aryl  C-H  oscillators,  this  is  equation  (15),  and  (16)  for  the 
methyl  C-H  oscillators, l5,|g 

NHH,(0l2/^l^H|0>2  (15)  |//J  =  (yU0l2(0L^I^2»93)|0)||0)2|0)3  (16) 


where  the  state  with  the  quantum  number  u  indicates  the  oscillator  that  absorbs  the  IR  radiation, 
and  the  o~0  states  indicate  ground  state  wave  functions.  The  relative  line  intensity  is  then  given 
by  equation  (17)  and  is  combined  with  the  Lorentzian  line  shape  function  to  simulate  the  1R 
spectra. 

Intensity  oc  (oyjn  ( 1 7) 
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w  here  (Do  is  the  transition  frequency,  and  T  is  the  full  width  half  max  (fwhm)  of  the  absorption 
feature.  Figure  4  shows  the  simulated  spectrum  of  the  C-H  stretches  in  the  TNT  molecule  for  a 
80  cm  line  width,  compared  with  the  near  1R  spectrum.  Fundamental  and  the  first  overtone 
transitions  are  simulated  w  ith  the  HCAO  model,  for  single  TNT  molecules  and  clusters  of  TNT. 
Figure  1A  shows  the  simplest  cluster  model  simulated  for  a  single  eight  molecule  unit  cell,  and 
Figure  1  B  shows  the  geometry  for  an  eight  unit  cell  cluster  of  TNT  molecules.  The  stick  models 
represent  the  molecules  modeled  with  UFF,  and  the  molecule  shown  with  the  ball  and  stick 
model  represents  the  molecule  that  was  treated  with  quantum  mechanics.  The  van  dcr  Waals  and 
electrostatic  forces  are  allowed  to  perturb  the  quantum  mechanical  model. 
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Figure  1 .  Cluster  models  for  neat  TNT.  (A)  is  a  single  molecule  surrounded  by  a  single  unit 
cell,  and  (B)  is  surrounded  by  eight  unit  cells.  The  ball  and  stick  model  represents  the  molecule 
treated  with  quantum  meehanies,  and  the  stick  models  represent  the  molecules  treated  with  the 
classical  Universal  Force  Field. 

2.2  Laboratory  Measurements 

Laboratory  measurements  provide  validation  of  the  computational  model  of  the 
molecule-surface  interaction.  The  validation  is  provided  by  agreement  of  predicted  overtone 
absorption  frequencies  to  experimental  measurements.  The  direct  measurement  of  overtone 
absorption  is  accomplished  by  photoaeoustie  speetroseopy.  When  molecules  absorb  optical 
energy,  they  have  an  increase  in  internal  energy  and  are  out  of  equilibrium  with  the  environment. 
The  most  efficient  pathway  to  equilibrium  with  the  environment  is  by  emitting  the  excess  energy 
as  heat.  When  the  heat  is  released  into  fluidie  surroundings,  such  as  air  or  solvent,  a  pressure 
wave  is  created  that  can  be  detected  with  an  acoustie  transducer.  4  The  measurement  is  directly 
proportional  to  the  absorbed  energy  and  is  completely  independent  of  light  scattering  from  the 
sample,  so  that  sample  morphology  has  little  effect  on  the  measurement  of  the  absorption.  0 
When  the  light  source  is  modulated,  the  resulting  acoustic  response  will  have  the  same 
modulation  frequency  with  some  delay. 

One  problem  with  a  photoaeoustie  approach  that  must  be  overcome  is  that 
overtone  transitions  tend  to  have  cross  sections  one  to  two  orders  of  magnitude  less  intense  than 
the  fundamental  transitions.  Comparison  between  the  methods  of  fundamental  and  overtone 
absorption  show  that  this  is  surmountable.  Although  measurement  of  fundamental  transitions 
occur  in  the  mid  to  far  IR  portions  of  the  spectrum  where  laboratory  light  sources  are  generally 
weak,  useful  measurements  are  routinely  obtained  in  the  laboratory.  Overtone  transitions  tend  to 
oeeur  in  the  near  IR  and  visible  wavelengths,  which  contain  more  powerful  light  sources.  Weak 
signals  are  overcome  with  more  intense  light  sources.  In  systems  that  arc  easily  saturated  with 
optical  energy,  the  measurement  is  performed  with  additional  sample  to  provide  measurable 
photoaeoustie  signal.’1 
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The  laser  photoacoustic  spectra  of  50  mg  of  neat  TNT  are  recorded  using  the 
MTEC300,  a  commercial  photoaeoustie  cell.  Tunable  excitation  radiation  of  about  0.6  mJ/pulsc 
is  produced  by  the  idler  of  a  Spectra  Physics  10  Hz  GCR-250  Nd:YAG  pumped  OPO  laser 
system  (Spectra  Physics  MOPO  730),  and  introduced  through  a  window  in  the  cell.  The 
wavelength  is  scanned  from  435  nm  (5180  cm'1)  to  465  mn  (6663  cm’1).  Because  the  excitation 
source  is  a  pulsed  laser,  the  acoustic  signal  is  also  pulsed,  although  the  acoustic  signal  has  a 
~1  ms  width,  compared  to  the  2.5  ns  pulse  width  of  the  laser.  Signal  is  acquired  with  a  Stanford 
Research  boxcar  integrator  with  its  gate  width  maximized  at  3  ps,  to  maximize  signal  to  noise. 
The  output  of  the  boxcar  is  digitized  with  a  National  Instruments  PCI-4472  digitizer,  and  the  data 
is  acquired  with  Labview  8.1.  Because  the  signal  is  not  modulated,  a  resonance  condition  for 
signal  amplification  is  not  required  nor  is  additional  digital  signal  processing. 

Sample  characterization  is  accomplished  via  photoacoustic  FT1R  measurements 
with  the  same  cell  that  is  used  for  the  overtone  measurements.  For  reference,  neat  TNT  samples 
are  also  measured  with  a  Perkin-Elmer  near  lR-Vis  spectrometer,  and  with  a  Nicolet  attenuated 
transmission/reflectance  FT1R. 


3.  RESULTS  AND  DISCUSSION 

3.1  Selection  of  High  Level  Method 

The  first  task  for  the  computational  component  is  the  determination  of  the 
appropriate  level  of  theory  for  treating  the  TNT  molecule  and  calculating  the  overtone 
frequencies.  To  obtain  the  overtone  frequencies  of  the  TNT  molecule,  selection  of  the 
theoretical  approach  is  reduced  to  methods  that  have  analytic  second  derivatives  of  the  energy 
with  respect  to  nuclear  position,  which  are  Hartree-Fock  (HF),  Density  Functional  Theory' 
(DFT),  and  Moeller/Plesset  Theory  (MP2/MP3).  The  anharmonic  terms  for  the  potential  energy 
surface  are  calculated  by  taking  the  numerical  third  derivative  of  the  energy  with  respect  to 
nuclear  position  along  the  normal  mode  coordinate.  The  STO-3G,  3-2 1G,  6-31+G(d.p),  and  6- 
311+G(d,p)  basis  sets  were  evaluated  for  geometry  optimization  and  vibrational  frequency 
calculation  utilizing  HF,  MP2,  and  DFT,  where  the  B3LYP  hybrid  functional  was  used  for  DFT. 
The  results,  such  as  the  fundamental  vibrational  frequencies  and  optimized  geometric  structure, 
were  compared  to  TNT  FTIR  spectra  measured  in  the  laboratory  and  published  X-ray 
crystallographic  data,  respectively.  The  best  results  arc  obtained  utilizing  B3LYP  and  MP2 
methods  with  the  6-31+G(d,p)  basis  set,  where  calculated  bond  lengths  and  angles  compare  well 
with  the  published  crystal  data,  and  calculated  C-H  stretching  vibrational  frequencies  compare 
well  with  the  laboratory  FTIR  measurements  (within  20  cm’1.  Figure  2).  The  HF  level  of  theory, 
without  scaling  factors,  produces  differences  of  ~225  cm’1  between  the  predicted  aromatic  C-H 
stretching  frequencies,  and  the  highest  frequency  in  the  FTIR  spectrum  (which  we  assigned  to 
the  aromatic  C-H  stretching  frequency).  For  the  first  overtone  frequencies,  the  same  trend  was 
also  followed  where  the  higher  levels  of  theory  resulted  in  lower  predicted  frequencies 
(Table  1).  A  small  peak  is  visible  at  6036  cm  1  in  Figure  2,  which  is  close  to  the  overtone 
frequency  for  the  aromatic  asymmetric  stretch  predicted  at  6161  cm  .  Although  the  B3LYP  and 
MP2  methods  proved  to  have  equivalent  accuracy  in  comparison  to  published  data  geometric 
data,  the  B3LYP  method  is  the  preferred  method  due  to  lower  computational  expense. 


17 


500  1000  1500  2000  2500  3000  3500 

Frequency  (cm1) 


Figure  2.  Comparison  of  calculated  harmonic  1R  absorption  spectra  and  a  TNT  film  absorption 
spectrum  measured  with  FT1R.  T  he  upper  trace  uses  the  Hartree-Fock  Method,  and  the  middle 
trace  uses  Density  Functional  Theory  with  the  B3LYP  functional.  Both  methods  use  the  6- 
31+G(d.p)  basis  set,  which  has  some  diffuse  functions. 

3.2  Comparison  of  FT1R  and  Near  IR  Data  to  Simulated  Spectra 

The  first  phase  of  experimental  data  involved  collection  of  TNT  spectra  with 
commercial  analytical  equipment.  The  mid  IR  absorption  spectrum  of  TNT  was  recorded  with 
an  Attenuated  Transmittance  Reflectance  (ATR)  spectrometer  with  the  data  shown  in  Figure  3. 
This  allowed  the  spectrum  of  a  thin  film  of  TNT  to  be  recorded.  Although  several  peaks  are 
present  in  the  recorded  spectrum  that  are  not  visible  in  the  published  IR  spectrum,  which  appear 
to  indicate  the  presence  of  water  and  CCfi  in  the  sample  chamber,  there  is  excellent  agreement 
for  the  C-H  stretch  frequencies  around  3100  cm  .  It  w'as  also  attempted  to  view  the  overtone 
absorptions  using  a  Perkin  Elmer  near  IR  spectrometer.  The  resulting  spectrum  is  shown  in 
Figure  4.  Although  a  peak  shows  up  around  ~3 1 00  cm  in  this  spectrum,  corresponding  to  the 
aromatic  C-H  stretches,  the  signal  to  noise  starts  to  become  rather  poor  at  ~5000  cm  where  the 
light  source  becomes  weak.  There  appears  to  be  a  small  peak  or  divot  in  the  spectrum  occurring 
about  where  quantum  mechanical  calculations  predict  the  overtone  absorption.  This  may  be  an 
overtone  absorption  of  neat  TNT. 

Figure  3  compares  the  results  of  the  standard  frequency  analysis,  with  anharmonie 
corrections,  within  Gaussian  and  HCAO  spectral  simulation  compared  to  the  photoacoustic  FTIR 
measurement  of  the  same  TNT  sample  used  for  laser  photoacoustic  measurement  of  the  overtone 
absorption.  As  shown  in  the  selection  of  the  level  of  theory  and  basis  set,  the  HF  predicted 
spectrum  differs  the  most  from  the  FTIR  measurement.  The  aryl  hydrogen  stretching  component 
of  the  predicted  spectrum  from  B3LYP  is  even  closer  to  experiment,  yet  the  contribution  from 
the  methyl  hydrogens  is  quite  different  than  the  FTIR  measurement.  Surprisingly,  given  that  the 
primary  assumption  of  HCAO  is  that  the  local  mode  description  is  valid  only  for  overtones,  the 
best  fit  is  obtained  from  the  HCAO  predicted  spectrum  of  from  single  TNT  molecules.  The  last 
simulated  spectrum  is  from  a  TNT  molecule  within  a  64  member  TNT  cluster.  The  model 
predicts  that  the  contribution  from  the  methyl  C-H  stretches  becomes  the  dominant  spectral 
feature.  This  is  in  complete  disagreement  with  FTIR  and  near  IR  measurements.  The 
disagreement  between  experimental  spectra  and  cluster  predicted  spectra  is  even  stronger  for  the 
overtones. 


18 


2000  2200  2400  2600  2800  3000  3200  3400  3600  3800  4000 

frequency  fcm-1) 


Figure  3:  Comparison  of  calculated  1R  absorption  spectra  (normalized  to  maximum  intensity)  and  a 
solid  TNT  powdered  sample  measured  with  photoacoustic  FT1R.  (A)  and  (B)  are  spectra  predicted 
by  frequency  analysis  of  the  normal  modes  generated  with  the  mass  weighted  Cartesian  coordinates 
and  anharmonic  corrections  applied  at  the  HF/6-31+G(d,p)  and  B3LYP/6-3 1+G(d,p)  levels  of  theory 
(C)  and  (D)  are  spectra  predicted  w  ith  the  HCAO  model  for  a  single  TNT  molecule  and  a  TNT 
molecule  in  an  8  unit  cell  cluster  of  TNT  (64  molecules)  at  the  B3LYP/6-31+G(d,p)  level  of  theory. 
Trace  E  is  the  photoacoustic  FT1R  measurement. 


The  experimental  overtone  spectrum  of  TNT  is  compared  to  several  predicted 
spectra  using  the  HCAO  model.  The  HCAO  model  as  shown  in  the  literature  uses  the  Morse 
oscillator  wave  functions  as  the  basis  set  for  the  calculation  of  transition  intensity.  When 
potential  energy  surface  (PES)  calculations  are  performed  for  a  TNT  molecule  within  a  cluster  of 
TNT  molecules,  the  PES  for  the  C-H  stretches  becomes  more  harmonic.  Given  that  the  other 
molecules  within  the  cluster  of  the  model  confine  the  hydrogen  along  the  stretching  coordinate, 
the  result  is  that  the  depth  of  the  potential  well  is  greatly  increased.  The  increase  in  well  depth 
makes  the  PES  down  at  the  bottom  more  harmonic.  Because  of  this,  as  shown  in  Table  2,  the 
anharmonicitics  predicted  by  the  HCAO  model  become  small.  The  exponential  term  in  equation 
(13)  becomes  so  small  for  small  anharmonicity  values  that  the  values  cannot  he  represented  by 
double  precision  floating-point  numbers.  Thus,  it  is  not  possible  to  compensate  vs  ith 
normalization  of  the  function,  and  it  became  necessary  to  use  linear  combinations  of  harmonic 
oscillator  wave  functions.  The  B  spectrum  in  Figure  5  shows  the  predicted  spectrum  from  the 
harmonic  oscillator  basis  set.  Instead  of  an  exact  match  between  A  and  B,  the  harmonic 
oscillator  wave  functions  lead  to  a  prediction  of  greater  intensity  for  the  methyl  C-H  stretches. 
This  is  also  in  disagreement  with  the  experimental  measurement,  which  does  not  show  a  strong 
methyl  C-H  contribution.  Contributions  from  higher  quantum  number  wave  functions  may  he 
necessary  to  bring  about  agreement  from  the  two  methods.  Although  the  overtone  spectrum 
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from  a  single  TNT  moleeule  model  with  the  Morse  wavefunction  basis  set  appears  to  agree  with 
the  experimental  spectrum  quite  well,  agreement  between  predieted  speetra  from  TNT  within 
clusters  and  the  experimental  data  beeomes  worse.  As  the  eluster  beeomes  bigger,  the  methyl  C- 
H  stretches  become  increasingly  dominant.  Correct  representation  of  the  wave  functions  may 
resolve  this  issue. 

There  arc  additional  considerations  that  may  account  for  increasing  discrepancy 
between  predicted  spectra  from  dusters  and  experimental  measurements  as  shown  in  Figure  5. 
The  cluster  models  so  far  have  only  examined  the  spectral  contribution  from  single  molecules, 
instead  of  an  average  over  an  ensemble  of  molecules.  Thus,  the  predicted  spectra  may  be  overly 
sensitive  to  the  local  environment  of  a  single  moleeule.  This  is  especially  true  of  the  predicted 
intensities,  or  transition  dipoles,  as  shown  in  Tabic  3,  where  the  aryl  C-H  transitions  are 
predicted  to  have  much  less  intensity  relative  to  the  methyl  C-H  stretches.  The  averaging  of 
several  predieted  moleeular  spectra  may  be  neeessary  to  adequately  predict  speetra.  An 
additional  consideration  is  that  the  van  der  Waals  and  electrostatic  interaction  terms  in  the  UFF 
force  field  may  cause  on  overestimation  of  how  much  the  surrounding  molecules  in  the  cluster 
perturbs  the  moleeule  under  study.  The  non-bonding  interactions  between  atoms  are  considered 
pairwise.  As  a  result,  three  body  effects,  such  as  the  screening  of  one  atom  by  another  are  not 
accounted  for.  It  may  be  necessary  to  introduce  sharper  cutoff  distances  for  the  non-bonding 
contributions  to  the  force  field,  or  to  increase  the  size  of  the  quantum  meehanical  model  to 
include  surrounding  molecules.  The  hybrid  moleeular  mechanics/quantum  mechanics  treatment 
also  uses  the  qEq  algorithm  to  assign  partial  charges  to  the  UFF  layer  in  the  model.  It  is  known 
that  this  algorithm  has  a  margin  of  error  for  this,  and  may  be  over  predieting  the  effeet  of  eharge 
on  the  TNT  molecule  in  the  quantum  layer  of  the  calculation.  Expansion  of  the  quantum 
mechanical  layer  will  allow  the  partial  eharges  on  the  surrounding  moleeules  to  be  determined 
directly  from  ab  initio  methods. 
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Figure  4:  Calculated  fundamental  and  first  overtones  of  the  C-H  stretching  absorptions  for  TNT 
using  the  HCAO  model  vs.  near  1R  transmission  measurement  of  a  TNT  film.  The  inset  shows  a 
small  feature  in  the  near  1R  spectrum  that  appears  at  -6000  cm'1,  on  a  predicted  overtone  feature. 
The  predicted  overtone  intensity  relative  to  the  fundamental  is  about  5%. 

With  respect  to  the  standard  anharmonic  frequency  analysis  within  the  Gaussian 
code,  except  for  very  small  clusters,  calculations  involving  standard  anharmonic  analysis  do  not 
successfully  complete.  Considerable  effort  was  devoted  to  determining  the  cause  of  this  failure. 
Although  the  geometry  optimization  of  the  model  system  typically  completes,  problems  arc 
encountered  during  the  frequency  analysis.  These  calculations  were  rerun  at  the  High 
Performance  Computing  Center  at  Aberdeen  in  an  attempt  to  diagnose  the  issue,  using  -100 
cores  and  the  LAN  parallel  version  of  Gaussian  2009.  A  similar  failure  occurred,  where  the 
calculation  did  not  complete  after  several  weeks  of  execution.  We  suspect  that  the  core  issue  is 
the  method  Gaussian  uses  to  determine  the  anharmonic  corrections  to  the  normal  mode 
frequencies.  To  determine  the  anharmonic  corrections,  Gaussian  needs  to  find  the  analytic 
representation  of  the  second  order  mass  weighted  force  constants  in  Cartesian  coordinates,  to 
find  the  third  order  force  constants.  For  models  involving  large  numbers  of  nonbonding 
interactions,  the  necessity  of  analytic  second  order  force  constants  may  be  an  impossible 
constraint.  In  contrast,  within  the  HCAO  model,  determination  of  the  third  order  force  constant 
merely  involves  fitting  a  fourth  order  polynomial  to  a  PES  along  a  single  coordinate.  As  a  result. 
HCAO  analysis  is  performed  on  systems  that  could  not  be  vicw'cd  with  the  standard  analysis. 
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Table  1 .  Predicted  harmonic  and  anharmonic  vibrational  frequencies  for  C-H  stretches  in  isolated 
and  clustered  TNT  molecules.  Columns  arc  in  order  of  decreasing  frequencies,  where  columns  1 
and  2  involve  the  aromatic  hydrogens,  and  3-5  correspond  to  motions  of  the  methyl  hydrogens. 
HCAO  Morse  and  HCAO  Harmonic  are  HCAO  models  using  either  Morse  oscillator  wave 
functions  or  harmonic  oscillator  wave  functions  as  basis  sets.  Linear  combinations  of  harmonic 
oscillators  are  necessary’  for  the  HCAO  analysis  when  the  anharmonicity  constants  become  small. 


C-H  Stretching  Frequency(cnf') 


Single  Molecule 

l 

2 

3 

4 

5 

HF/6- 3 1  +  G (cl p)  A tt harmonic 

3321.21 

3314.09 

3205.23 

3202.74 

31  11.26 

1st  Overtone 

6580.28 

6570.02 

6349.26 

6349.26 

6362.47 

B3L  Y P/6-3 l+G(d,p)  Anharmonic 

31 14.60 

3110,40 

3003.40 

2963.00 

2963.00 

1st  Overtone 

6172.00 

6161.50 

5977.10 

5960.60 

5878.00 

B3L  Y P/6- 31 + G(d,p)  HCAO  Morse 

3061.50 

3061.50 

2855.20 

2887.20 

2887.53 

1st  Overtone 

6039.39 

6039.39 

5660.11 

5727.34 

5728.06 

B3 L  YP/6- 3  J  +G(d p)  H CA O  harmonic 

3061.50 

3061.50 

2855.20 

2887.20 

2887.50 

1st  Overtone 

6039.40 

6039.40 

5660.10 

5727.30 

5728.10 

Molecular  Cluster 

HF/6- 3 1  +  G(d,p)  Anharmonic  1  unit  cell 

3325.14 

3318.83 

3312.34 

3207.11 

3150.32 

1st  Overtone 

6550.39 

6538.50 

6530.42 

6371.49 

6286.35 

B3LYP/6-31+G(d,p)  1  unit  cell  Anharmonic 

3124.68 

3120.81 

3116.98 

3014.98 

2987.63 

1st  Overtone 

6138.36 

6125.90 

6123.00 

5959.52 

5954.17 

B3L  YP/6- 3 1 + G(d.p)  l  unit  cell  HCAO 

3006.35 

3002.73 

2901.39 

2812.44 

281 1  47 

1st  Overtone 

6012.71 

6005.45 

5842.06 

5661.82 

5659.94 

B3L  YP/6-3J +G(d.p)  8  unit  cell  HCAO 

3011.00 

3011.40 

2812.90 

2817.90 

2801.20 

1st  Overtone 

6022.00 

6021.10 

5658.80 

5673.20 

5636.30 
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Tabic  2.  Morse  parameters  for  HCAO  models  of  an  isolated  TNT  molecule,  a  TNT  molecule 
within  a  single  unit  cell  (eight  molecules  total),  and  a  TNT  molecule  within  an  eight  unit  cell 
cluster  of  TNT  molecules  (64  molecules  total).  Note  that  for  the  models  involving  clusters,  the 
anharmonieity  constant  becomes  small.  This  effect  affects  the  method  of  calculating  the 
absorption  intensities. 


Morse  Parameters 

I 

B3I }  P/6-3  1 +G(d.p )  (ot  (cm  1 ) 

3229.20 

Anharmonieity  (o^\c  (cm  ’ ) 

83.73 

B3L  VP  6-3  f+G(d.p)  /  unit  cell  (ot  (cm  ' ) 

3006.40 

Anharmonieity  tol.\r  (cm  ) 

0.01 

B3L )  P  6-3  /  +  G(dp)  8  mut  celt  a\  (cm  1 ) 

3010.30 

Anharmonieity  (cm'1) 

0.23 

2 

3 

4 

5 

3229.20 

3094.50 

3131.10 

3131.30 

83.73 

87.72 

87.45 

87.41 

3002.80 

2971.70 

2893.20 

2892.10 

0.02 

0.87 

4.77 

4.72 

3010.40 

2893.80 

2896.70 

2887.00 

0.09 

6.65 

4.00 

7.54 

Table  3.  Calculated  transition  dipoles  in  units  of  Debye/eni.  The  first  two  rows  are  for  single 
TNT  molecules  using  the  Morse  and  harmonic  wave  function  basis  sets,  respectively.  Models 
with  a  large  molecular  cluster  around  the  TNT  molecule  used  only  harmonic  oscillator  wave 
function  basis  sets  to  calculate  intensity. 


12  3  4  5 


B3LYP  6-31  +  Gulp)  HCAO  Morse 

2.954713 

2.954539 

0.259592 

0  663128 

0.664081 

1st  Overtones 

0.143255 

0.143281 

0.193207 

0  046579 

0.046634 

B3L  YP  6-3  /  +  G(dp)  HCA  O  HB 

132.7 

132.7 

34.039 

22  293 

22.281 

1st  Overtones 

1.6258 

1 .6257 

0.14619 

0.14635 

0.14622 

B3L  YP  6-3 l+GUlp)  l  cell  HCAO 

3.501575 

3.617291 

0. 1 83668 

1  719555 

3.03683 

1st  Overtoms 

0.010915 

0.010844 

0.15329 

0.074628 

0.066534 

B3LYP  6-3  f +G(d.p)  8  cell  HCAO 

7.7497 

2.9994 

0.50615 

0.010906 

0.000812 

1st  Overtones 

0  001 149 

0.000108 

1.1818 

1.4X92 

0.07721 
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4. 


CONCLUSIONS 


Spectroscopic  measurements  of  the  fundamental  and  first  overtone  C-H  stretching 
are  measured  with  photo  acoustic  FTIR,  photoacoustic  laser  spectroscopy,  transmission  FTIR 
spectroscopy,  and  near  1R  transmission  spectroscopy.  Single  TNT  molecules  as  well  as  several 
TNT  eluster  systems  are  investigated  with  theoretical  chemistry  with  the  objective  generating 
predicted  IR  spectra  from  the  C-H  stretching  modes  for  comparison  to  experiment.  HCAO  and 
standard  frequency  analysis  arc  applied  to  models  of  single  TNT  molecules,  eight  molecule 
clusters  of  TNT,  and  64  moleeule  dusters  of  TNT.  Only  the  HCAO  analysis  successfully 
predicted  overtone  spectra  for  comparison  to  experimental  measurement.  The  standard 
anharmonie  frequency  analysis  has  been  unsuccessful  in  predicting  speetra,  most  likely  due  to 
limitations  with  the  core  method  of  determining  the  anharmonie  corrections  themselves. 
Although  HCAO  analysis  has  successfully  produced  predicted  speetra,  the  spectra  for  clusters 
are  not  yet  in  good  agreement  with  spectroscopic  measurements.  The  disagreement  may  arise 
either  from  ( 1 )  an  overestimation  of  the  nonbonding  interactions  between  molecules  in  the 
cluster;  (2)  fact  that  the  spectra  are  calculated  from  single  molecules  rather  than  an  ensemble  of 
molecules;  (3)  because  the  boundary  between  the  quantum  mechanical  and  the  classical 
mechanical  models  has  not  been  well  chosen;  or  (4)  because  an  inadequate  basis  set  of  harmonic 
oscillator  wave  functions  has  been  used  in  the  intensity  calculations. 

The  issues  identified  above  will  be  resolved  as  follows.  When  using  the 
nonbonding  interactions  of  a  force  field  to  perturb  the  quantum  mechanical  portion  of  the  model, 
the  contribution  from  each  atom  in  the  molecular  mechanics  layer  is  considered  pairwise.  The 
effect  of  this  approach  is  to  ignore  screening  effects  from  other  atoms,  and  overestimate  the  total 
nonbonding  interaction.  Additional  calculations  will  be  run  that  introduce  steeper  cutoffs  for 
nonbonding  interactions  in  the  force  field  model.  Another  possible  source  of  inaccuracy  may  be 
poor  selection  of  the  boundary  between  the  quantum  mechanical  layer  and  the  molecular 
meehanies  model.  This  ean  be  remedied  with  additional  calculations  that  surround  the  target 
molecule  w  ith  an  additional  layer  of  quantum  mechanical  molecules.  We  expect  more  accurate, 
partial  charges  on  the  neighboring  molecules.  Another  consideration  is  the  contribution  of 
multiple  moleeules  to  the  overtone  spectra.  As  stated  previously,  the  overtone  spectrum  is 
sensitive  to  the  local  environment,  and  as  a  result,  each  molecule  that  contributes  to  the  speetrum 
may  have  a  mueh  different  individual  spectrum.  To  address  this  issue,  additional  calculations 
will  explore  how  the  speetrum  ehanges  if  it  results  from  an  average  over  several  molecules  in  the 
model.  F  inally,  if  the  above  approaches  do  not  resolve  discrepancies  between  predieted  overtone 
speetra  from  dusters,  the  effeet  of  the  harmonie  oseillator  basis  set  in  the  intensity  calculation 
will  be  investigated.  Higher  quantum  number  wave  functions  will  be  included  in  the  linear 
combination. 


The  stated  objective  of  this  study  is  to  1)  validate  computational  models  by 
agreement  between  experimental  and  predicted  overtone  spectra,  and  2)  to  apply  the  model  to  the 
test  case  of  TNT  on  alumina.  To  obtain  these  objectives,  additional  work  needs  to  be  done. 
First,  additional  overtone  spectra  of  TNT  on  alumina  need  to  be  measured,  where  the  primary 
spectral  contribution  is  from  TNT  molecules  interacting  with  the  surface.  Second,  additional 
modifications  to  the  HCAO  model  are  necessary  to  produce  better  agreement  between  predicted 
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speetra  and  experiment.  Lastly,  onee  predieted  speetra  have  better  agreement  with  experiment, 
to  use  the  models  to  prediet  interaction  properties,  sueh  as  binding  energy. 


Simulated  Spectra  r 


Figure  5.  Comparison  of  laser  PAS  speetaim  of  TNT  compared  to  HCAO  simulated  spectra. 
Trace  (A)  is  the  speetrum  obtained  from  a  single  TNT  moleeule  and  a  Morse  Oscillator  basis 
set,  (B)  is  the  same  molecule  with  a  Harmonie  oseillator  basis  set,  (C)  is  a  T  NT  molecule  in  a 
1  unit  cell  eluster,  and  (D)  is  a  single  moleeule  within  an  8  unit  eell  cluster,  where  all  predicted 
speetra  have  been  normalized  with  respect  to  the  most  intense  peak. 
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